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Abstract 

We study the behavior of the optimal path between two sites separated by a distance r on a 
(i-dimensional lattice of linear size L with weight assigned to each site. We focus on the strong 
disorder limit, i.e., when the weight of a single site dominates the sum of the weights along each 
path. We calculate the probability distribution P(£ opt \r, L) of the optimal path length £ opt , and 
find for r < L a power law decay with £ opt , characterized by exponent g opt . We determine the 
scaling form of -P(^ pt| r > L) in two- and three-dimensional lattices. To test the conjecture that the 
optimal paths in strong disorder and flow in percolation clusters belong to the same universality 
class, we study the tracer path length of tracers inside percolation through their probability 
distribution P(£ tT \r, L). We find that, because the optimal path is not constrained to belong to a 
percolation cluster, the two problems are different. However, by constraining the optimal paths 
to remain inside the percolation clusters in analogy to tracers in percolation, the two problems 
exhibit similar scaling properties. 
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I. INTRODUCTION 



Flow in porous media, due to its ubiquitous nature, has received a great deal of attention 
in recent decades . Interest has been driven by multiple real-world problems 

_water polution 0, Q]. Porous media is typically 
and flow through a variety of models ranging 



such as oil extraction j^, |5( and grounc 
simulated through percolation systems 



rom t ic idealized an1 in a. labyrinth to the very practical ones like eonveetiye tracer Mow 



mi. 

In the study of flow in percolation clusters, it has been hypohesize that transport proper- 
ties may be related to some geometric property of the cluster. However, despite considerable 
effort [J such a geometric property relation has not been found jg]]. Here, we argue that, 
for the particular case of convective tracer flow in percolation clusters, a connection exists 
between transport properties and the static properties of a model of optimal paths in dis- 
ordered lattices, described in detail below. This relation allows for a mapping between the 
two problems that opens new approaches for their study. 

The first hint of a possible connection between flow in porous media and the optimal path 
problem arouse in the work of Lee et al. j3| . This reference deals with a simplified model of 
flow in porous media related to secondary oil extraction, in which an invading fluid (water, 
steam, etc.) injected at A pushes oil out through extraction well B. This method is used 
when the oil reservoir does not have enough pressure to be exploited without this added 
driving pressure. The system is considered in the following way: the medium is modeled 
by a percolacion cluster at the critical threshold density, and a steady state laminar flow 
of an incompresible fluid is established between A and B. Then, tracers which mimick the 
flow of the driving fluid are injected at A and travel purely by convection to B (Sec. IIII|) . 
The tracer length (£ tr ) probability density function (PDF), P(£ tr ), was measured and it was 
found that the most probable traveling length £l r , defined as the value of £ tr for which P(£tr) 
is maximum, scales with the distance r between the injection and extraction sites A and B 
as 

£* r ~r d * r [dtr = 1.21 ± 0.02 for d=2]. (1) 

The remarkable feature of Eq. (0) is that, for d = 2, the exponent d tr « d opt , where d opt is 
the exponent for the optimal path length in the strong disorder limit as defined by Cieplak 
et al. I3|(see below). In this study, we propose that transport in percolation can be 
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directly related to the optimal path in strong disorder, and we support our proposal with 
extensive numerical simulations for <i-dimensional lattices with d = 2 and 3. 



The optimal path problem formulated by Cieplak in Ref. 13] is that of finding the path 
of lowest cost to go from one end of a (i-dimensional lattice to the other end when to each 
site (or bond) i, we associate a weight e, = e aXl , with x t G [0,1). This is equivalent to 
choosing 6j from the distribution 

W{e) = — ee[l,expo]. (2) 
as 

The energy of any path of length £ on the lattice is given by the sum 

i 

E = J2 e ^ ( 3 ) 

3 

where j is an index running over the sites of the path. The limit a — > oo is known as the 
strong disorder limit. The optimal path of length £ opt is the path for which E is minimal 
with respect to all other paths. The optimal path length scales with r as [k| 



? opt ~ r*** 



[dopt = 1-22 ±0.01]. (4) 



The first reference that we have found to this problem, although formulated in a different 
context is that of Ambegaokar Q], where the metal-insulator transition is considered by 
a percolation model in which sites represent possible electron states which can be reached 
only by hopping through quantum tunneling, and the hopping rates are exponential, mo- 
tivating Eq. (J2J). Another context in which these weights have been observed is that of 
magnetoresistance in thin Ni films (if] ]. 

Our work can be equated to the following question: is flow in percolation clusters, a 
dynamical process, connected to the optimal path length in strong disorder, a static property, 
as suggested by the similarity of Eqs. (JTJ) and (jU? To this end, we study the PDF P(£ opt \r, L) 
for the optimal path to have a length £ opt , given a system size L, and an Euclidean distance 
r between the starting and ending sites A and B of the path. We will compare this PDF 
with P(£ tT \r, L), the PDF that convective tracer paths have a length £ tr in a percolation 
system of size L at criticality, where the starting and ending sites are at a distance r. 

There are some indications that this connection is indeed present, given that other re- 
lations between percolation and strong disorder optimal paths have been reported. For 
instance, for a lattice with disorder given by Eq. O) in the strong disorder limit, the most 



probable largest weight of the site used by the optimal path is e apc , where p c is the perco- 



lation threshold [14|. Also, Wu et al. 17j recently determined through the study of tracer 
flow on a lattice with disorder, that the strong disorder limit has a length scale that scales 
as a u , where v is the connectedness exponent of percolation [3], and hence a system is in 
the strong disorder limit only when a u > L. In this paper, we study the ultrametric limit 
a — > oo. 

In Sec. HU we present results for the distribution P(£ opt \r, L). In Sec. IHIl results for 
P(£ tT \r,L) are presented and compared with P(£ opt \r, L). We then discuss the results in 
Sec. ED 



II. OPTIMAL PATH DISTRIBUTION 

To study P(£ opt \r, L) we use the "bombing algorithm" proposed in Ref. The optimal 
path length £ opt between sites A and B is found by eliminating (bombing) sites of the lattice 
in decreasing order of weight, but leaving those sites necessary to keep A and B connected. 
When all sites that do not disconnect A and B are eliminated, only the sites of the optimal 
path remain 3- 

In Fig. Elwe present P(£ opt \r, L) for r < L for a square lattice of sites. Four distinct 
features appear: 



The most probable optimal path length £* pt scales with r as 



£* . ~ r d ° 

'-opt ' 



(5) 



The values of d opt have been reported for several lattice dimensions d, and also have 



been shown to be universal 



2l\ . Here, we rescale P(£ opt \r, L) with the exponent d opt , 



calculated elsewhere for the average optimal path length £ opt , but, as Figs. [T^i and[TjD 
show, d opt also produces the correct scaling for £* t . Our results for d opt are reported 
in Table H] for d = 2 and 3. 



A lower cutof ? 
percolation 
form 



Fig . H T b)1 which, in analogy with the distribution of minimal paths in 
22j, |23j], is expected to be a stretched exponential function f 1 of the 



exp(— ax 



x = 



^-opt 



(6) 
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where a is a lattice-dependent constant, and opt is a universal exponent satisfying 

0,0 



"opt 1 



An upper cutoff due to the effect of the finite lattice size L. A stretched exponential 
behavior is also expected to describe this region J, Q through a function f 2 of 
the form 



Mv) = exp(-Py*°*) 



y 



*-opt 

L d °p t 



where (3 is a lattice-dependent constant, and ijj opt has universal propert 



A power-law region described by 



ies J. 



19 



(8) 



231. 
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P{L»t) ~ Ct° pt t rd ° pt < %t < (9) 
The above considerations lead us to postulate for P(£ opt \r, L) a full scaling Ansatz 
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P{£ op t\r,L) 



-opt 



-ffopt 



/l 



^•opt 

where the prefactor l/r^p* is necessary for normalization. We have tested this Ansatz for 
d = 2, 3 and found it to be consistent with our earlier simulations jsij ] . 

An interesting feature of P(£ opt |r, L) is that, as d increases, g opt decreases. In other 
words, the longer optimal paths at larger dimensions have a larger probability (see Table [IJ). 
Additionally, since g opt < 2 for all d, £ op t and all higher moments diverge as L — > oo. 

To calculate the exponents opt and ip opt of Eqs. (jOJ and (JBJ), we introduce the function 



h 



"opt 



(10) 



n 



^opt 



,X,A 



In 



.4 



P(£ op t|r,L)r rf °p 
which, upon using Eqs. (jBJ), (JHJ) and ffTOfl yields 



3o P t 



(11) 



n 



--opt 



/7 t ' ' 



A, A 



In 



.4 



, f/opt\ f (to 



In A + a 



-opt 



\-rfopt 



'°Pt \-dopt 



V'opt 



(12) 



We have made use of A = L/r in the argument of the function f 2 so that M^opt/L dopt ) = 
M ^ ~ dopt ^opt/ r dopt )• The constant A is an auxiliary parameter chosen to make the minimum 



value of II slightly larger than unity. Defining x = £ op t/r dopt , we show in Fig. El n(x, A, A) 
for d = 2 and the fit lines for the exponents of both fx and fi-, which are reported in Table UJ 
The values of opt we calculate are close to the values predicted by Eq. (|ZJ) for d = 2, 3. 



III. COMPARISON BETWEEN FLOW IN PERCOLATION AND THE OPTIMAL 
PATH 

We now study the PDF P(£ tI \r, L) with the purpose of comparing it to P(£ op t\ r , L), and 
address analyze the detailed conditions under which optimization and flow in percolation 
occur. Our analysis (see Sec. IIV|) explains the differences we observe, and also the "right 
way" in which the two problems become equivalent. 

Since simulations for flow on percolation clusters are performed, we describe the two- 
dimensional case of the algorithm [5]. We represent the reservoir as a two-dimensional 
site percolation cluster, and choose sites at (— r/2,0) and (r/2,0), denoted by A and B, 
respectively, to be the injection and extraction well positions. Points A and B are separated 
by a geometric distance r, and the system box has corners at (±L/2, ±L/2). We construct 



24|. 



percolation clusters at p c using the Leath algorithm 

To model tracer motion we use the analogy with electrical circuits, where for each bond, 
the pressure drop corresponds to the voltage difference, and the flow corresponds to the 
electrical current on the bond. A pressure difference between sites A and B drives the 
tracer. For each realization, 10 4 tracers are introduced at site A, and then collected at site 
B. The set of all sites through which there is a non-zero current defines the cluster backbone 
of M B sites. 

The "pressure" difference across bonds is equivalent to a "voltage" difference, so by 
solving Kirchhoff 's equations on the backbone, we obtain the potential (pressure) drops AV 
over all bonds for a given realization. Due to the nature of the model, no turbulence or 
other complex fluid flow effects are considered, which is equivalent to assuming laminar flow 
inside the system. Additionally, given that the bonds have vanishing radius, the tracer flow 
"perfectly mixes" at the nodes. For site i having Sj outgoing bonds, the tracer selects a 
bond with a probability 

w H = Y~k\- b' = !>•••>*; i = l,...,M fl ]. (13) 
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For incoming bonds, Wij = 0. This guaranties that the tracer dynamics are completely 



convective, i.e., with infinite Peclet number 



26| |. The total traveling length of a tracer 



is the number of bonds of the path connecting A and B, chosen by this tracer. Since the 
particles do not interact with one another, it is equivalent to launching one particle at a 
time into the cluster. This procedure is known as particle launching algorithm jloL fll | . We 
determine the probability distribution of the tracer traveling lengths P(£ tv \r, L) by counting 
the number of particles that travel from site A to site B along a path of length £ tT , over all 
the particles and all realizations of the percolation cluster. 

The form of P(£% T \r y L) for the two-dimensional case was suggested in Q]. Here, we 
extend these results to d = 3 (Fig. [HJ) • Once again, the power law and stretched exponential 
behaviors are present. The scaling of the most probable tracer path length is given by 



These results yield 



tr 



tr 



Ldt: 



(14) 



u 



<-tr 



Hr 



L du 



(15) 
(16) 



where functions h\ and hi have the forms 

hi(z) = exp(-/x,z" 0tr ) 

and 

h,2(u) = exp(—pu^ tT ) 

Arguments similar to those leading to Eq. (fT2*j) indicate how to determine the exponents </> tr 
and iptr, reported in Table |TJ 

The power law region is characterized by the exponent g tr , which is different from g opt 
(see Table H}. We present in Fig. 0] curves for the power law regime for both P(£ tr \r, L) and 
P(£ opt \r, L) in d = 2 and 3. In Table H] we see the difference in the slope of the power law 
decay between P(£ tT \r,L) and P(£ opt \r, L). Moreover, as d increases, g opt decreases and g tT 
increases, indicating the differing behaviors for the two problems. In the next section, we 
explain the origin of the differences, how these differences can be removed, and under which 
conditions the two problems coincide. 



IV. DISCUSSION 



The numerical results presented above show the difference in the values of the scaling 
exponents g opt and g tr of the distributions. To understand these differences, we now elaborate 
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on the characteristics of the optimal path problem in comparison to those of tracer paths in 
percolation. 

In Fig. EJa) we represent the optimal path in strong disorder, where the dark areas 
represent regions with site weights e« = e aXi with x-i < p c , and the white areas regions 
with site weigths e = e aXi with Xj > p c . Typically, the arbitrary choice of A and B may 
lead to a path connecting them that requires visiting regions with site weights e > e apc . 
In contrast, the tracers inside percolation clusters must, by definition, travel on the same 
percolation cluster (spanning or otherwise), because the flow takes place only if there is a 
percolating path between A and B. Therefore, this difference between the flow and optimal 
path problems presents a possible explanation for the differences between P(£ opt \r, L) and 
P(£ tT \r, L). Optimal paths tend to be longer because they are able to visit more sites of the 
lattice and are therefore of longer length, whereas tracers in percolation flow are confined to 
a given cluster, and their traveling lengths are much more limited. These features intuitively 
explain why g tr is larger than g opt . 

The above considerations lead to the following hypothesis: if the optimal path search is 
constrained to pairs of sites within regions of the lattice that are part of the same cluster 
[Fig. 03(b)], then the scaling of P(£ opt \r, L) and P(£ tT \r, L) would coincide. To test this, we 
present P(£ opt \r, L) and P(£ tr \r,L) in Fig. where the optimal paths satisfy the condition 
that their highest weight is at or below percolation. This condition forces the optimal paths 
to be inside percolation clusters. Indeed, for this case (Fig. EJ), the two quantities exhibit 
very similar behavior, supporting our hypothesis. The exponent g op t inside percolation 
now becomes very close to g tr . On the other hand, d op t does not change, confirming the 
equivalence of the two problems. We also have similar results for three-dimensional lattices. 

In summary, we have shown that P(£ opt \r, L) has a power law tail with an exponent 
g opt which decays as d grows and is different from the power law tail of P(£ tr \r, L). This 
difference seems to be related to the fact that the optimal path crosses percolation clusters 
and thus tends to have longer lengths compared with tracers which are always inside perco- 
lation clusters. When £ opt is measured only inside percolation clusters, our results suggest 
P(£ op t\ r , L) and P(£ tr \r, L) are equivalent and the two problems possibly belong to the same 
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universality class. 
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TABLE I: Exponents characterizing P(£ op t\r,L) and P(ltr\ r , L), which are defined in the text. 
The value of g op t is determined from P(£ pt\r = 4, L = 256) in Fig. ^a), for which the power law 
region is the longest. The values of 4> op t and ^opt are from Fig. 2. 



Optimal Path in Strong Disorder ("Static") 


d 


"opt 


<?opt 


<^opt (calculated) 




V'opt (calculated) 


2 


1.22 ± 0.01 [13., 21] 


1.55 ±0.05 


4.8 ±0.5 


4.55 ±0.21 


5.3 ±0.3 


3 


1.42 ± 0.02[27J 


1.37 ±0.05 


2.1 ±0.1 


2.3 ±0.1 


4.3 ±0.3 


Optimal Path in Strong Disorder inside Percolation ("Modified Static") 


2 


1.21 ±0.02 


1.82 ±0.05 


4.9 ±0.4 


4.76 ± 0.45 


2.3 ±0.4 


3 


1.40 ±0.03 


2.2 ±0.1 


2.0 ±0.1 


2.5 ±0.2 


3.6 ±0.2 


Tracer Path ("Dynamic") 


d 




9tx 


0tr (calculated) 


^ = dtr _i 


■(/V (calculated) 


2 


1.21 ± 0.02[12j 


1.82 ±0.05 


4.7 ±0.4 


4.76 ± 0.45 


2.7 ±0.2 


3 


1.37 ±0.05 


2.23 ±0.09 


1.81 ±0.02 


2.7 ±0.4 


3.46 ± 0.04 
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1 10 



Scaled optimal path length l Jr° vt 

FIG. 1: (a) Distribution P(£ opt \r, L) for r = 4 and system sizes L = 64, 128 and 256. As L increases, 
the power law region with exponent g opt becomes better defined, and the upper cutoff increases, 
(b) Probability distribution P(£ opt \r,L) for (r = 8,L = 64) (solid line), (r = 16, L = 128) (dashed 
line), and (r = 32, L = 256) (dotted line) for two-dimensional systems. The ratio between L and 
r is kept fixed for these curves, (c) Scaled distribution r dopt P(£ op t\r, L) vs. scaled optimal path 
length t ov t/f dopt for the curves in (b). The collapse has been achieved using the exponent d opt 
reported for £ opt , which is also valid for the most probable length i* opt as evidenced in the plot. 
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L=64, r=8 
128,16 
256,32 
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Scaled optimal length x=l jr opt 



FIG. 2: The scaling function U(x, A, A) for ^4 = 0.1 vs. the scaled optimal path length x = l pt/r dopt 
for system sizes (r = 8, L = 64), (r = 16, L = 128), and (r = 32, L = 256). The two straight lines 
serve as guides to the eye for the data that determine the exponents 4> opt and ip op t- 



13 




FIG. 3: Probability distribution P(£ tT \r, L) for r = 2 and L = 128 in d = 2 and d = 3. In a similar 
fashion as for P(£ opt \r, L) we see a power law region that we characterize by exponent gtr- Another 
feature of this plot is the increasing steepness of P(£ tr \r, L) as d increases (implying g tT increases 
with d), a feature for which P(£ opt \r, L) has the opposite behavior, as g op t decreases with d. 
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FIG. 4: The power law tails of P(£ tl \r,L) (solid) and P(£ opt \r,L) (dashed) in d = 2 and 3. The 
upper pair is for d = 2 with L = 256 and r = 4. The lower pair is for d = 3 with L = 128 and 
r = 2. 
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opt 



FIG. 5: (a) Schematic of occupied sites for p below the percolation threshold p c , and the optimal 
path in strong disorder. The darker regions represent sites that are still present when p c is reached. 
We see in this case that £ op t must cross the region above p c (i.e., leave the cluster) to connect A 
and B. (b) If sites A and B are chosen within the a cluster below p c , the optimal path does not 
leave the cluster because that would increase the cost. 
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FIG. 6: Comparison of P(£ op t\r,L) inside percolation with P(£ tr \r,L) for (r = 16, L = 64), (r = 
32, L = 128) and (r = 64, L = 256). The solid lines represent the optimal path distribution, and 
the long dashed lines the tracer length distributions. The values of r and L have a fixed ratio equal 
to L/r = 4. The similarity between distributions is clear, supporting our hypothesis. The small 
separation along the horizontal axis between P(£ op t\r, L) and P(£t T \r, L) (consistent for the three 
pairs of curves) is due to non-universal details of the two models. 
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